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A  new  method  is  given  for  preventing  the  simplex  method  from  cycling. 
Key  features  are  that  a  positive  step  is  taken  at  every  iteration,  and  nonbasic 
variables  are  allowed  to  be  slightly  infeasible.  There  is  no  additional  work  per 
iteration.  Computational  results  are  given  for  the  first  53  test  problems  in 
netlib,  indicating  reliable  performance  in  all  cases. 

The  method  may  be  applied  to  active-set  methods  for  solving  nonlinear 


programs  with  linear  constraints. 
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1.  Introduction 


Degeneracy  is  often  regarded  as  a  discomforting  but  otherwise  tolerable  hindrance 
to  the  simplex  method,  and  to  other  active-set  algorithms  for  solving  optimization 
problems  involving  linear  constraints.  Sequences  of  non-improving  steps  are  known 
to  occur  (perhaps  many  times  during  a  run),  but  such  sequences  are  rarely  observed 
to  be  infinite.  The  phenomenon  of  “stalling”  is  therefore  recognized  and  accepted, 
but  “cycling”  is  deemed  very  unlikely  to  occur. 

In  spite  of  such  folklore,  a  rigorous  anti-cycling  procedure  can  provide  welcome 
peace  of  mind  to  users  and  implementors  alike,  particularly  if  the  cost  is  small. 
Such  a  procedure  was  given  by  Wolfe  [Wol63],  and  the  possible  benefits  have  been 
demonstrated  recently  by  Ryan  and  Osborne  [R086].  (See  also  Falkner  and  Ryan 
[FR87].)  Our  aim  here  is  to  present  an  alternative  anti-cycling  procedure  that,  like 
Wolfe’s  method,  involves  little  overhead  and  has  proved  to  be  effective  in  practice. 
We  also  investigate  the  relationship  with  Wolfe’s  method. 

*The  material  contained  in  this  report  is  based  upon  research  supported  by  the  Air  Force  Office 
of  Scientific  Research  Grant  87-01962;  the  U.S.  Department  of  Energy  Grant  DE-FG03-87ER25030; 
National  Science  Foundation  Grants  CCR-8413211  and  ECS-8715153;  and  the  Office  of  Naval  Re¬ 
search  Contract  N00014-87-K-0142. 
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An  anti-cycling  procedure 


1.1.  The  standard  LP  problem 


Most  of  our  discussion  will  be  in  terms  of  the  simplex  method  [Dan63]  and  the 
primal  linear  programming  problem 


minimize  cTx 

subject  to  Ax  —  b,  l  <  x  <  u, 


(1) 


where  A  £  9?mxn  (to  <  n).1  Let  Ax  =  Bxb  +  Nx /v  denote  the  partitioning  of  A 
and  x  into  basic  and  nonbasic  variables,  where  B  £  9?mX7n  js  the  usual  basis  matrix. 

In  typical  .iTnplcmcr.tatior.s  ^thc  simplex  method,  the  general  constraint"  Ax  —  b 
are  satisfied  throughout  and  infeasibility  occurs  only  with  respect  to  the  bounds. 
For  most  iterations,  factors  of  the  current  B  are  obtained  by  updating,  but  at  the 
beginning  of  a  run  and  periodically  thereafter,  B  is  factorized  directly  and  the  basic 
variables  are  recomputed  to  satisfy  Bxb  +  Nxpi  =  b.  If  the  newly  computed  xb  does 
not  satisfy  its  upper  and  lower  bounds  to  within  some  feasibiliiy  tolerance  6  >  0, 
Phase  1  of  the  simplex  method  is  invoked  to  move  the  infeasible  variables  towards 
their  violated  bounds.  Phase  2  starts  (or  resumes)  when  the  feasibility  tolerance  is 
satisfied. 

A  similar  optimality  tolerance  is  used  to  judge  whether  any  reduced  costs  are 
sufficiently  positive  or  negative  to  give  an  improved  solution.  Note  that  the  feasibil¬ 
ity  and  optimality  tolerances  are  typically  of  order  10-6,  which  is  much  larger  than 
a  typical  machine  precision  e  («  10-16). 

In  practice  this  means  that  “optimal”  solutions  are  in  reality  feasible  and  near- 
optimal  solutions  to  the  perturbed  problem 

minimize  cTx 

subject  to  Ax  =  6,  fg  —  be  <  .Tg  <  ttg  +  6e,  (2) 

In  <  xn  <  w-as 


where  e  is  a  vector  of  ones,  and  B  and  N  relate  to  the  final  basis  obtained. 

For  the  anti-cycling  procedure  developed  here,  we  find  it  convenient  to  allow 
nonbasic  variables  to  be  infeasible  in  a  similar  way.  Thus,  n  place  of  (1)  and  (2), 
we  aim  to  find  a  feasible  and  near-optimal  solution  to  the  perturbed  problem 

minimize  < 'Jx 

(3) 

subject  to  Ax  =  b,  l  -  be  <  x  <  u  +  be. 


Since  the  final  B-N  partition  is  somewhat  unpredictable,  we  anticipate  that  prac¬ 
titioners  accustomed  to  (2)  will  find  problem  (3)  essentially  equivalent.  In  practice, 
few  (if  any)  nonbasic  variables  will  terminate  infeasible. 

In  terms  of  conventional  error  analysis,  the  constraints  of  (3)  can  always  be 
satisfied  numerically,  as  long  as  6  is  sufficiently  large  compared  to  machine  precision. 

'Implicitly,  A  and  x  are  of  the  form  A  =  (  A  I ).  x  =  ( i  s),  where  s  is  a  set  of  slack  variables 
with  appropriate  bounds  in  l  and  u.  However,  we  never  need  to  distinguish  between  £  and  s. 
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(This  is  true  even  if  (1)  has  no  feasible  solution.)  Note  however  that  we  don’t 
rigorously  “use  up”  all  the  freedom  allowed  by  S  to  optimize  the  objective  function; 
i.e.,  if  a  solution  appears  to  be  optimal  we  don’t  try  to  improve  the  objective  at  the 
expense  of  making  additional  variables  slightly  infeasible.  Instead,  we  are  content 
to  terminate  knowing  that  some  basic  or  nonbasic  variables  may  lie  up  to  6  outside 
their  true  bounds. 

1.2.  Key  aspects 

Two  main  features  are  involved  in  the  new  anti-cycling  procedure: 

•  The  feasibility  tolerance  S  is  increased  slightly  at  the  start  of  every  iteration. 

•  Numerical  values  are  stored  for  all  components  of  x. 

The  first  feature  allows  a  positive  step  to  be  taken  at  every  iteration.  The  second 
allows  slight  infeasibilities  to  be  recorded  correctly  when  basic  variables  become 
nonbasic.  We  shall  speak  of  the  “EXPAND”  procedure2  when  referring  to  this 
particular  refinement  of  the  simplex  method. 

The  EXPAND  procedure  is  a  practical  row-selection  method,  in  the  sense  that 
it  chooses  a  pivot  row  in  the  simplex  method.  We  are  able  to  retain  the  “maximum 
pivot”  property  of  the  row-selection  method  due  to  Harris  [Har73].  In  addition, 
we  are  able  to  remove  an  unforeseen  weakness  in  previous  implementations  of  the 
Harris  procedure. 

1.3.  Other  anti-cycling  methods 

Wolfe’s  method  [Wol63]  and  the  lexicographic  rule  [DOW55,Dan63]  are  both  row- 
selection  procedures.  As  with  the  EXPAND  procedure,  these  methods  can  be  used 
with  any  column-selection  (“pricing”)  strategy.  In  particular,  they  allow  “partial 
pricing” — an  important  advantage  when  n  >  m. 

Recently,  a  new  column- selection  method  has  been  given  by  Dantzig  [Dan88]  and 
further  developed  by  Klotz  [Klo88].  It  may  be  used  with  any  row-selection  method. 
An  additional  “pricing”  vector  is  required  and  the  method  is  not  directly  amenable  to 
partial  pricing.  However,  it  appears  to  be  promising  for  highly  degenerate  problems. 

The  first  anti-cycling  method  of  Bland  [Bla77]  prescribes  both  the  pivot  column 
and  the  pivot  row,  so  again  does  not  allow  partial  pricing. 

In  contrast,  Benichou  et  al.  [BGHR77]  perturb  the  vector  b  and  then  apply  a 
“normal”  simplex  procedure.  If  the  perturbation  is  chosen  randomly,  the  probability 
of  cycling  is  zero.  The  EXPAND  procedure  is  somewhat  akin  to  this  approach, 
particularly  when  the  perturbation  is  removed  (Section  4.4).  However,  we  do  not 
rely  on  random  numbers,  and  we  recommend  a  considerably  smaller  perturbation 
than  the  O(10~3)  referred  to  in  [BGHR77], 

The  primal-dual  methods  of  Balinski  and  Gomory  [BG63],  Graves  [Gra65j  and 
Fletcher  [Fle85]  involve  a  nested  sequence  of  subproblems  similar  to  those  arising  in 
Wolfe’s  method. 

2  EXPanding-tolerance  ANti-Degeneracy  procedure 
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An  anti-cycling  procedure 


2.  Nonbasic  Solutions 

Historically,  the  simplex  method  has  been  implemented  in  such  a  way  that  numerical 
values  are  recorded  for  the  basic  variables  only.  The  value  of  each  nonbasic  variable 
is  typically  implied  (by  a  status  indicator)  to  be  the  variable’s  lower  or  upper  bound. 
Often  the  bounds  are  implicitly  zero  and  infinity  (0  <  Xj  <  oo  for  all  j)  and  nonbasic 
variables  are  implicitly  zero. 

In  more  general  implementations,  some  or  all  variables  Xj  are  allowed  to  have 
arbitrary  bounds  lj  and  uj,  and  the  value  of  a  nonbasic  variable  is  normally  defined 
to  be  one  of  the  bounds;  thus,  Xj  =  lj  or  Uj  according  to  a  status  indicator.  A 
complication  arises  with  “free  variables”  (satisfying  -oo  <  x y  <  +oo),  since  provi¬ 
sion  must  be  made  for  them  to  be  nonbasic  even  though  they  are  likely  to  be  basic 
at  a  solution.  Typically,  a  nonbasic  free  variable  is  defined  to  have  the  value  uro, 
thereby  avoiding  the  need  to  store  any  otner  value.  (This  approach  was  used  in 
various  versions  of  MINOS  up  to  and  including  MINOS  5.0  [MS83] .) 

History  aside,  many  implementors  have  recognized  that  certain  benefits  arise 
if  an  arbitrary  value  can  be  stored  for  each  and  every  variable.  For  example,  in 
mathematical  programming  systems  such  as  MPSX/370  and  MPS  III,  the  BASIC 
procedure  is  designed  to  input  numerical  values  for  any  number  of  variables  and 
produce  from  them  a  basic  solution.  Further  examples  are  the  “pseudo-constraints” 
of  Fletcher  and  Jackson  [FJ74],  the  “temporary  constraints”  of  Gill  and  Murray 
[GM78],  and  the  “pegged  variables”  of  Nazareth  [Naz86,Naz87],  The  essential  idea 
is  that  nonbasic  variables  can  be  temporarily  frozen  at  specified  values. 

Thus  in  linear  and  nonlinear  programming,  while  the  term  nonbasic  is  often 
taken  to  mean  “equal  to  zero”  or  “equal  to  an  upper  or  lower  bound”,  it  is  more 
useful  to  define  a  nonbasic  variable  as  one  that  is  currently  fixed  at  a  specified  value. 
The  working-set  strategy  involved  will  eventually  allow  such  a  variable  to  move  as 
if  it  were  being  released  from  a  normal  bound. 

This  definition  was  adopted  in  MINOS  5.1  [MS87].  Explicit  values  are  stored  for 
all  variables,  and  at  each  iteration  of  the  simplex  and  reduced-gradient  algorithms, 
nonbasic  variables  are  allowed  to  take  any  strictly  feasible  value: 

lj  <  Xj  <  uj.  (4) 

An  advantage  during  cold  starts  is  that  variables  can  be  initialized  at  the  “safe” 
value  of  zero  in  cases  where  a  user  has  specified  deceptively  large  bounds,  such  as 
lj  =  —10s,  it y  =  +108.  (There  is  no  need  to  initialize  ij  at  one  of  its  bounds.) 
Similar  advantages  arise  when  restarting  modified  problems  and  recovering  from 
singular  bases. 

For  the  anti-cycling  procedure  of  this  paper,  we  have  generalized  (4)  by  allowing 
nonbasic  variables  to  be  slightly  outside  their  bounds: 

lj  —  6  <  Xj  <  Uj  +  6,  (5) 

where  6  is  the  feasibility  tolerance  mentioned  earlier.  This  also  eliminates  a  difficulty 
with  the  Harris-type  steplength  procedure,  as  we  now  describe. 
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3.  Steplength  Procedures 

In  general,  optimization  algorithms  proceed  by  generating  a  search  direction  p  and 
then  changing  the  variables  according  to  x  <—  x  +  op  for  some  steplength  a  (a  >  0). 

In  the  primal  simplex  method,  x  always  satisfies  the  constraints  Ax  =  b.  The 
“pricing”  or  “column- selection”  strategy  chooses  a  nonbasic  variable  to  be  moved 
from  its  current  value.  (This  variable  usually  enters  the  basis.)  A  search  direction 
p  is  then  determined,  such  that  A(x  +  op)  =  b  for  any  step  a. 

The  subsequent  steplength  computation  is  known  as  the  “ratio  test”  or  the  “row- 
selection”  procedure.  The  aim  is  to  find  which  variable  is  the  first  to  encounter  a 
bound.  (This  variable  usually  leaves  the  basis.) 

3.1.  The  standard  ratio  test 

A  “textbook”  ratio  test  assumes  that  the  current  point  x  is  feasible  (/  <  x  <  it)  and 
finds  the  largest  step  a  that  keeps  the  new  point  feasible:  l  <  x  +  ap  <  u.  Some 
blocking  variable  xr  reaches  one  of  its  bounds  exactly,  so  that  xr  +  apr  =  lT  or  ur 
depending  on  the  sign  of  pr.  For  each  j,  let  ay  be  the  step  that  takes  xy  to  one  of 
its  bounds: 

f  (h  -  *j)/pj  Pi  <  o» 

=  S  («>  -  Xj)/Pj  Pj  >  0.  (6) 

[  oo  otherwise. 

Further,  let  ay  =  miny  ay.  The  maximum  feasible  step  is  then  a  =  ar  >  0,  and  the 
blocking  variable  is  ir. 

Given  i,  p,  /  and  u,  we  see  that  the  ratio  test  determines  a  steplength  a  and  an 
index  r.  We  shall  refer  to  this  procedure  by  writing 

(a,r)  =  ratio-test  (x,p,  l,  u). 

A  danger  with  the  textbook  ratio  test  is  that  the  pivot  element  pr  could  be 
very  small  if  xT  is  close  to  its  bound.  (If  pr/||p|i  is  small,  the  next  basis  matrix 
will  be  ill-conditioned.)  To  provide  a  nominal  safeguard,  we  define  a  “cut-off”  value 
below  which  small  elements  py  are  treated  as  zero  in  (6).  (In  MINOS,  “small”  means 
|py|  <  tolp  where  tolp  =  e1^3  for  linear  programs  and  e2|/3||p||  for  nonlinear  programs, 
with  e2'3  «  10~n.  Some  implications  are  discussed  in  Section  5.1.) 

3.2.  The  Harris  ratio  test 

In  [IIar73],  Harris  observed  that  some  freedom  in  choosing  r  can  often  be  introduced 
by  using  a  two-pass  procedure.  The  first  pass  determines  a  perturbed  steplength  a  1 
that  is  slightly  too  large  to  keep  x  feasible: 

(al,rl)  =  ratio-test (x,p, l  -  6e,u  +  be ), 

where  6  is  the  usual  feasibility  tolerance  («  10-6).  It  is  important  to  note  that  x  is 
assumed  to  satisfy  the  perturbed  bounds  (/  -  Se,u  +  be).3  It  follows  that  al  >  0. 

3In  Phase  1  of  the  simplex  method,  some  components  of  l  and  u  are  altered  to  ±oo  (perhaps 
implicitly)  to  make  this  true;  see  Section  7. 
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An  anti-cycling  procedure 


Figure  1:  Paths  followed  with  standard  (a-b-c-d)  and  Harris  (a-f-d)  ratio  tests. 


The  second  pass  then  considers  all  unperturbed  steps  aj  (6)  no  larger  than  al, 
choosing  the  index  associated  with  the  largest  pivot: 

|pr|  =  maxjpy|  such  that  aj<al. 

We  define  the  corresponding  steplength  to  be  a2  =  ar .  The  step  a  =  a2  is  then 
acceptable  as  long  as  it  is  positive.  In  general,  we  define  a  =  max{a2,0}. 


3.3.  An  example  with  a2  positive 

Let  (/,u)  =  (0,oo)  and  suppose  the  current  feasible  solution  is  x  —  (0.0009, 1)T. 
The  Harris  ratio  test  with  6  —  10"3  would  then  lead  to  a  step  x  <—  x  +  ap  of  the 
form 

(  -0.0001  \  /  0.0009  \  (  -0.1  \ 

l  0  H  1  -looj- 

The  steplength  is  positive  {a  —  a2  —  0.01,  r  =  2).  This  is  slightly  too  large  to  keep 
x  feasible,  but  the  larger  pivot  is  successfully  chosen. 

Figure  1  illustrates  a  similar  case.  The  standard  ratio  test  would  lead  to  the  path 
(a-b-c-d);  the  solution  would  stay  strictly  feasible,  but  the  constraint  encountered 
at  b  corresponds  to  a  rather  small  pivot  element.  By  allowing  this  constraint  to  be 
slightly  violated,  the  Harris  ratio  test  would  choose  a  less  oblique  constraint  to  be 
encountered  at  f,  giving  the  shorter  and  numerically  more  reliable  path  (a-f-d). 
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3.4.  An  unexpected  error  < 

« 

Note  that  a2  will  be  negative  if  the  blocking  variable  lies  slightly  outside  its  bound  ' 

and  p  is  leading  it  further  from  the  same  bound.  (For  example,  xi  above  is  now  > 

—0.0001  and  p\  could  be  negative  at  the  next  iteration.)  . 

The  natural  inclination  is  to  set  a  =  0  and  interpret  this  as  a  zero  or  “degen-  1 

erate”  step.  Note  however  that  the  blocking  variable  becomes  nonbasic.  In  most  [ 

implementations  of  the  simplex  method,  this  means  that  the  blocking  variable  is  j 

actually  changed  (perhaps  implicitly)  to  lie  exactly  on  its  bound.  ; 

In  effect,  if  the  blocking  variable  is  slightly  infeasible  (say  lT  -  8  <  xr  <  lT),  then  j 

most  implementations  of  the  Harris  procedure  change  the  variables  according  to  i 

x  <-  x  +  aer,  (7) 

where  |a|  <  S.  Such  a  change  produces  an  unintentional  error  in  satisfying  Ax  =  b. 

The  error  may  be  of  order  8,  which  is  typically  much  larger  than  e. 

In  practice,  errors  of  this  kind  tend  to  be  eliminated  each  time  the  basis  is 
refactorized,  since  the  basic  variables  are  typically  recomputed  in  order  to  satisfy 
Ax  =  b  accurately.4  Provision  is  made  to  return  to  Phase  1  if  the  recomputed 
variables  lie  outside  their  bounds  by  more  than  8.  On  well-behaved  problems,  few 


iterations  (if  any)  are  required  to  regain  feasibility,  but  in  runs  lasting  thousands  t 

of  iterations,  the  risk  of  a  few  extra  iterations  every  50  (a  typical  factorization  t 

frequency)  amounts  to  a  nontrivial  overhead.  In  the  worst  case,  “few”  can  be  more 
than  50  and  an  optimum  may  not  be  achieved. 

For  nonlinear  problems,  the  perturbation  to  x  in  (7)  can  cause  a  discontinuity 
in  the  nonlinear  functions  and  may  lead  to  a  failure  in  the  linesearch  procedure. 


3.5.  A  simple  cure 

To  avoid  this  difficulty  with  the  Harris  ratio  test,  it  would  be  sufficient  to  implement 
a  “zero”  step  literally.  If  a  blocking  variable  is  slightly  infeasible,  we  should  make 
it  nonbasic  and  retain  its  infeasible  value  rather  than  moving  it  onto  its  bound.  The 
variable  should  be  temporarily  frozen  at  that  value  (across  basis  factorizations  if 
necessary)  until  the  normal  pricing  strategy  allows  it  to  move.  Provision  should  still 
be  made  to  revert  to  Phase  1  after  refaciorization,  but  given  a  .table  b»  sis-  hand  ling 
package,  the  likelihood  of  losing  feasibility  will  be  greatly  reduced. 

An  alternative  is  to  allow  a  negative  step  whenever  a2  <  0,  giving  the  blocking 
variable  a  chance  to  move  exactly  onto  its  bound.  This  approach  has  been  used  in 
the  quadratic  programming  and  linear  least-squares  codes  QPSOL  3.2  and  LSSOL 
1.0  ( GMSW84,GHM*86 ].  However,  it  is  then  necessary  to  perform  a  ratio  test  on 
the  reverse  search  direction  —  p,  obtaining  a  possibly  different  blocking  variable  that 
again  may  be  unable  to  reach  its  bound  exactly.  Since  the  objective  value  will  move 
slightly  in  the  wrong  direction,  there  is  also  the  possibility  of  cycling. 

We  propose  a  further  alternative  next. 


*Solve  Bib  =  b  —  Nxn,  or  preferably  solve  By  =  b  —  Ax  and  update  ib  •—  ip  +  y 
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4.  An  Anti-cycling  Procedure 


One  way  to  prevent  cycling  of  the  classical  kind  is  to  ensure  that  zero  or  negative 
steps  never  occur.  In  the  proposed  procedure  we  insist  that  a  >  0,  so  that  the 
objective  function  always  improves. 

Given  a  feasibility  tolerance  6,  suppose  as  before  that  the  current  x  is  feasible 
to  within  that  tolerance: 

/  -  6e  <  x  <  u  -f  6e.  (8) 

Now  suppose  that  6  is  changed  to  a  slightly  larger  tolerance  8  as  follows: 


6  =  6  +  r,  where  0  <  r  <  <5. 


(9) 


Since  6  >  8,  we  have 


l  -  8e  <  x  <  u  +  8e, 


(10) 


and  it  is  clearly  possible  to  take  a  positive  step  in  any  direction  p  before  encountering 
a  perturbed  bound.  To  find  an  acceptable  positive  step,  we  apply  the  Harris  ratio 
test  in  the  normal  way.  If  the  resulting  step  a2  is  negative  or  too  small,  we  replace 
it  by  a  certain  step  a^n  as  follows. 


EXPAND  procedure: 

1.  (First  pass)  Define  a  “largest  allowable  step”  al  using  the  increased  feasibility 
tolerance  8: 

( al,rl )  =  ratio-test (x,p, l  -  6e,u-\-  8e). 

2.  (Second  pass)  Find  the  step  a2  and  index  r  associated  with  the  largest  allow¬ 
able  pivot: 

a2  =  ar  where  \pr  j  =  max  \pj\  such  that  a3<al. 

3 

(The  quantities  a}  are  defined  by  equation  (6).) 

3.  Define  a  “smallest  allowable  step”  am;n  =  T/\Pr\- 

4.  If  a2  >  amin,  set  a  =  o?.  (When  x  changes  to  x  +  ap,  this  step  allows  the 
blocking  variable  xT  to  reach  its  bound  exactly.) 

5.  Otherwise,  set  a  =  Omjn.  (In  this  case,  the  new  value  of  xT  will  “overshoot” 
its  bound,  but  its  infeasibility  will  be  no  greater  than  8.) 

The  first  pass  gives  a  positive  step  al,  and  from  (8)-(  10)  it  is  clear  that 

al  >  t j\prl\. 

Since  the  second  pass  maximizes  the  pivot  element,  we  then  have 
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and  it  follows  that  amin  is  both  positive  and  not  too  large.  It  is  preferable  to  take 
tne  step  a2  whenever  possible  (to  allow  xT  to  reach  its  bound),  but  if  a 2  is  negative 
or  too  small,  the  step  Omit,  is  acceptable. 

To  summarize,  we  define  a  =  rnax{a2,amin}  and  x  =  x  +  ap.  We  have  shown 
that  a  positive  step  is  taken  (am;n  <  a  <  al)  and  that  the  new  point  satisfies  the 
required  bounds: 

l  —  be  <  x  <  u  +  6.  (11) 

Since  (11)  is  analogous  to  (8),  the  process  can  be  repeated  once  the  feasibility 
tolerance  is  (again)  increased  as  in  (9). 

4.1.  Infeasible  nonbasic  values 

Let  A  be  the  distance  between  the  blocking  variable  and  the  corresponding  bound: 
A  =  \xr  -  lT |  or  |xr  -  ur|. 

When  A  >  r,  the  preferred  step  a  =  a2  is  taken  and  the  blocking  variable 
reaches  its  bound  exactly.  This  is  the  normal  “nondegenerate”  case. 

If  A  <  r,  the  step  a  =  om;n  is  taken  and  the  blocking  variable  moves  a  total 
distance  of  r  and  terminates  infeasible.  The  final  value  of  xr  must  be  recorded  when 
the  blocking  variable  is  made  nonbasic. 

Figure  2  illustrates  a  normal  step  and  two  examples  of  a  degenerate  step.  We 
assume  the  blocking  variable  xr  is  being  constrained  by  its  lower  bound  (so  pr  <  0). 
The  sloping  lines  plot  the  value  of  xr  +  apT  against  a,  with  three  possible  starting 
values  for  xr.  The  lower  bound  is  at  the  horizontal  a  axis. 

In  the  top  case,  xT  is  relatively  large  initially  and  reaches  its  bound  after  a 
reasonably  large  step.  We  take  the  normal  Harris  step  a  =  a2  >  amjn. 

In  the  middle  case,  xr  starts  out  feasible  but  reaches  its  bound  after  a  very  small 
step.  We  insist  on  taking  a  larger  step  a  =  Qmin,  and  the  variable  becomes  slightly 
infeasible.  We  count  this  as  a  “degenerate”  step,  even  though  a  positive  move  is 
made. 

In  the  third  case,  xr  is  already  infeasible  and  becomes  even  more  infeasible. 
However,  after  a  step  o  =  omin  it  still  satisfies  the  required  bound  xr  >  lT  —  6. 
Again  we  count  this  as  a  degenerate  step. 

The  interpretation  of  “a  degenerate  step”  is  that  “a  slight  infeasibility  has  just 
been  created  among  the  nonbasic  variables.  The  objective  function  has  improved 
in  compensation”.  Since  it  is  common  for  blocking  variables  to  reenter  the  basis 
at  some  later  iteration,  the  total  number  of  nonbasic  infeasibilities  at  any  stage  is 
generally  less  than  the  number  of  degenerate  steps  so  far. 

4.2.  Typical  parameter  values 

The  preceding  sections  have  discussed  the  steplength  computation  for  one  iteration 
of  the  simplex  method.  Various  parameters  are  involved  in  a  complete  implementa¬ 
tion,  as  listed  below.  We  indicate  specific  values  that  might  be  used  in  practice  on 
a  machine  with  about  16  digits  of  precision. 
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XT  +  apr 


Figure  2:  Three  possible  starting  values  for  the  blocking  variable  xr. 

6  =  10~6  is  the  “main”  feasibility  tolerance. 

So  —  0.55  is  the  feasibility  tolerance  used  at  the  start  of  a  cycle  of  iterations. 

=  0.995  is  the  feasibility  tolerance  reached  at  the  end  of  a  cycle  of  iterations. 

K  —  10000  is  the  number  of  iterations  in  a  cycle. 

r  =  (6x  -  6o)/K  «  1 10“ 10  is  the  amount  by  which  the  feasibility  tolerance  is 
incremented  each  iteration. 

6k  =  6, t_i  +  r  is  the  feasibility  tolerance  used  at  the  fc-th  iteration  of  a  cycle  (k  =  1 
to  K). 

In  the  present  implementation,  the  “main”  tolerance  6  may  be  set  by  default  or 
specified  by  the  user.  Note  that  6q  and  6k  are  both  similar  to  5.  The  aim  is  to 
keep  the  “current”  feasibility  tolerance  5*  much  the  same  as  the  one  intended  by  the 
user,  but  to  increase  it  steadily  from  5o  to  6k  over  a  rather  long  cycle  of  consecutive 
iterations. 
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Figure  3:  Expansion  of  feasible  region  as  6k  increases  from  to  Sr. 

4.3.  Illustration 

Figure  3  depicts  a  feasible  region  that  expands  between  the  two  dashed  lines  as  <S k 
increases  from  f0  to  6jc  (K  =  5).  The  first  feasible  point  is  at  vertex  a,  and  the  step 
towards  vertex  b  will  not  go  beyond  the  first  dotted  boundary. 

If  the  feasibility  tolerance  were  zero,  the  simplex  method  would  follow  the  path 
(a-b-c-d-e-f).  With  a  positive  tolerance,  the  step  (c-d)  would  be  lengthened  and 
a  slight  infeasibility  would  arise  temporarily,  as  illustrated  in  Figure  1.  Vertices  b, 
c  and  f  would  be  reached  exactly. 

The  path  indicated  by  arrows  might  be  taken  if  the  feasible  region  were  defined 
by  certain  additional  hyperplanes  (not  shown).  These  hyperplanes  would  be  in 
higher  dimensions  and  would  have  to  cause  near-degeneracy  at  each  of  the  expanded 
vertices  corresponding  to  b,  c,  d  and  e.  The  main  idea  is  that,  even  if  an  iterate  lies 
on  or  close  to  a  degenerate  vertex  of  the  current  feasible  region  (i.e.,  close  to  one 
of  the  dotted  lines),  a  forward  move  will  always  be  possible  within  the  next  feasible 
region  (defined  by  the  next  dotted  line). 

4.4.  A  resetting  procedure 

At  certain  stages  we  require  a  “resetting  procedure”  to  remove  nonbasic  infeasibili¬ 
ties.  The  main  steps  are  as  follows. 

1.  The  values  of  nonbasic  variables  are  scanned.  Any  that  lie  within  b  of  a  bound 
are  moved  exactly  onto  the  bound.  (This  will  include  variables  that  were 
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slightly  infeasible  when  they  were  last  removed  from  the  basis.)  A  count  is 
kept  of  the  number  of  nontrivial  adjustments  made — say,  those  greater  than 
10"10. 

2.  If  the  count  is  positive,  the  basic  variables  are  recomputed  from  the  remaining 
variables,  thereby  satisfying  Ax  =  6  to  (essentially)  machine  precision. 

3.  The  current  feasibility  tolerance  is  reinitialized  to  So¬ 
li  a  problem  requiies  more  than  K  iterations,  we  invoke  the  resetting  procedure 

and  continue  with  a  new  cycle  of  K  iterations.  (The  decision  to  resume  in  Phase  1 
or  Phase  2  is  based  on  $j.) 

We  also  invoke  the  resetting  procedure  when  the  optimizer  reaches  an  apparently 
optimal,  infeasible  or  unbounded  solution,  unless  this  situation  has  already  occurred 
R  times,  where  R  is  a  further  parameter.  If  any  nontrivial  adjustments  are  made, 
iterations  are  continued.  Typically,  R  -  1  would  be  sensible  for  apparent  optimality, 
since  in  most  practical  cases  the  optimality  test  is  satisfied  (again)  immediately  after 
the  reset.  The  final  solution  is  then  “conventional”  in  the  sense  that  no  nonbasic 
variables  lie  outside  a  bound.  In  badly  conditioned  cases,  an  arbitrary  number  of 
iterations  may  be  needed  to  regain  feasibility  and  optimality  following  a  reset,  and 
R-  2  may  be  preferable,  since  the  second  reset  will  normally  adjust  fewer  nonbasics 
and  there  is  still  a  chance  of  terminating  at  a  “conventional”  solution. 

Note  that  the  solution  obtained  after  k  iterations  in  a  cycle,  or  k  iterations  after 
a  reset,  is  feasible  to  within  6k  (assuming  Phase  1  has  terminated).  We  may  regard 
all  preceding  iterations  as  a  means  of  reaching  such  a  point,  and  it  is  irrelevant  that 
the  feasibility  tolerance  has  been  adjusted  during  the  process.  Since  6k  <  6,  it  should 
be  acceptable  to  terminate  at  such  a  point  if  the  optimality  test  is  satisfied.  We 
can  therefore  advocate  using  R  =  0  if  there  is  concern  over  the  arbitrary  number  of 
iterations  that  may  be  required  following  a  reset.  In  other  words,  though  we  must 
reset  every  K  iterations,  there  is  no  real  need  to  reset  once  the  optimality  criteria 
are  satisfied.  The  solution  will  satisfy  the  constraints  of  problem  (3)  with  feasibility 
tolerance  equal  to  the  current  6k- 

A  similar  situation  exists  in  the  degeneracy-resolving  procedure  of  Benichou 
et  al.  [BGHR.77,  pp.  292-294],  in  which  a  perturbation  of  order  6  (6  =  10-2  or 
10-3)  is  added  to  the  right-hand  side  vector  b.  Once  the  perturbed  problem  has 
been  solved,  the  perturbation  is  removed  and  the  dual  simplex  algorithm  is  applied 
(often  requiring  no  further  iterations).  If  6  were  reasonably  small  (say  6  =  10-6),  one 
could  argue  that  the  solution  to  the  perturbed  problem  be  accepted  for  all  practical 
purposes. 

4.5.  Convergence 

In  ou.  use,  the  only  question  of  non- convergence  arises  with  resetting  every  K 
iterations.  If  K  were  small,  the  potential  loss  of  ground  after  each  reset  could 
conceivably  lead  to  a  classical  cycle  of  period  K.  Our  choice  of  K  —  10000  is 
intended  to  make  the  probability  of  such  a  cycle  negligible. 
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To  emphasize  the  point,  we  note  that  previous  implementations  of  the  sim¬ 
plex  method  have  been  operating  (in  effect)  with  K  set  to  the  basis  factorization 
frequency — typically  50  or  less.  Failures  due  to  cycling  have  been  rare  (though 
not  completely  absent;  for  example,  see  Benichou  et  a  1.  [BGHR77,  pp.  292-294]). 
Various  other  implementation  details  were  probably  contributing  factors. 

To  a  large  extent,  the  chance  of  failure  due  to  resetting  depends  on  cond(B), 
the  condition  number  of  a  typical  basis  matrix.  If  cond(B)  approaches  1/c  «  1016 
where  e  is  the  machine  precision,  then  any  algorithm  is  likely  to  fail.  However,  there 
should  be  no  risk  of  failure  when  cond(-B)  approaches  1/6  w  106.  By  choosing  K 
large  and  retaining  the  values  of  slightly  infeasible  nonbasic  variables  across  basis 
factorizations,  we  essentially  remove  all  risk. 

5.  A  Simplified  Procedure 

A  preliminary  implementation  of  the  EXPAND  procedure  was  used  for  the  experi¬ 
ments  conducted  by  Lustig  [Lus87].  This  version  was  simpler  and  potentially  more 
efficient  on  nondegenerate  problems;  we  therefore  describe  it  briefly.  As  before,  we 
assume  that  the  feasibility  tolerance  has  just  been  increased  to  6  =  S  +  r. 

Simplified  EXPAND  procedure: 

1.  (First  pass)  Obtain  (al,rl)  -  ratio-test(x,p,l,u). 

2.  Define  a  “smallest  allowable  step”  =  t/\Pti\. 

3.  If  al  >  Oniin,  set  (a,r)  =  ( al,rl )  and  exit. 

4.  (Second  pass)  Otherwise,  set  (a,r)  =  ratio-test  (x,p,  l  -  6e,u  +  6e). 

Ironically,  this  approach  reverses  the  two  passes  in  the  Harris  procedure.  It 
has  the  advantage  of  terminating  frequently  after  the  first  pass  (which  is  just  the 
classical  ratio  test  applied  to  the  true  problem  data).  The  blocking  variable  reaches 
its  bound  exactly. 

If  a  second  pass  is  required,  the  blocking  variable  must  be  made  nonbasic  at  a 
slightly  infeasible  value  ( lT  -6  or  ur  +  6). 

A  possible  disadvantage  is  that  the  pivot  element  \pT\  is  not  maximized  within 
a  set  of  candidates.  Nevertheless,  the  final  step  satisfies  a  >  r/|pr|  whether  one  or 
two  passes  are  required.  This  tends  to  prevent  selection  of  a  small  pivot  element, 
unless  the  feasible  region  is  unbounded.  We  do  not  expect  numerical  instability  if  S 
and  r  have  the  recommended  values  (Section  4.2).  No  difficulties  were  encountered 
in  the  computational  tests. 

5.1.  The  effect  of  ignoring  small  elements  of  p 

A  crucial  requirement  of  the  EXPAND  procedures  is  that  all  components  of  x  be 
at  least  a  distance  r  inside  the  current  perturbed  bounds  (1  —  6,  u  +  6).  If  small 
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elements  of  p  are  ignored  during  the  computation  of  a  (Section  3.1),  there  is  a  slight 
risk  that  the  required  property  will  not  hold  for  the  next  iteration. 

The  risk  exists  if  a\pj\  >  r  for  any  ignored  elements  pj ;  i.e.,  if  a  >  rftolp. 
For  typical  parameter  values,  this  means  if  a  >  5.  In  such  cases,  once  x  has  been 
updated  to  x  +  ap,  a  suitable  precaution  would  be  to  test  if  any  components  lie 
outside  the  bounds  (/  —  6,u  +  6).  Any  that  do  could  be  moved  onto  those  bounds. 

To  date  we  have  not  included  such  a  precaution  in  our  implementations.  It 
would  be  more  strongly  recommended  if  the  parameters  e,  6  and  r  were  substantially 
different  from  those  assumed  here. 

An  alternative  is  to  ignore  fewer  elements  of  p  (by  reducing  tolp),  since  when 
t  t  there  is  essentially  no  danger  of  a  small  pivot  being  selected  even  if  tolp  =  0. 
However,  it  is  common  for  many  elements  of  p  to  be  very  small,  and  excluding  such 
elements  from  the  ratio  test  can  give  significant  savings  on  large  problems. 

6.  Relationship  to  Wolfe’s  Procedure 

The  EXPAND  procedure,  simplified  or  otherwise,  may  be  interpreted  as  a  modifi¬ 
cation  of  Wolfe’s  anti-cycling  procedure  [Wol63],  as  we  will  now  show.  We  discuss 
the  case  with  general  lower  bounds  on  x  (but  no  upper  bounds). 

6.1.  Wolfe’s  procedure 

Let  LPo  denote  the  problem  to  be  solved: 

LPo  minimize  cTx 

subject  to  Ax  -  b,  x  >  /. 

Wolfe’s  “ad  hoc”  procedure  takes  effect  when  the  simplex  method  encounters  a 
degenerate  feasible  vertex  (say  x  =  xq).  The  degeneracy  structure  of  xo  is  used  to 
define  the  following  subsidiary  linear  program: 

LPi  minimize  cTx 

subject  to  Aa:  =  6,  xD  >  lD  —  d,  xN  >  lN, 

where  d  is  a  positive  vector,  D  denotes  basic  variables  that  are  currently  on  a  bound, 
and  N  denotes  variables  that  are  currently  nonbasic.  (Problem  LPi  is  the  same 
as  LPo  except  for  the  bounds  on  the  basic  variables.  For  degenerate  variables 
the  bounds  have  been  relaxed,  and  for  the  other  basic  variables  they  have  been 
removed.)5 

Clearly,  a:o  is  a  non-degenerate  feasible  point  for  the  new  problem.  When  the 
simplex  method  is  applied  to  LPj,  the  values  obtained  for  x  are  not  directly  relevant 
to  LPo,  but  the  bases  generated  (and  the  associated  dual  variables)  have  meaning 
for  both  problems.  Three  situations  may  arise: 

5Wolfe’s  procedure  can  be  described  using  alterations  to  x  and  6  as  well  as  l.  but  the  concepts 
are  the  same. 
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1.  A  finite  optimum  is  obtained  for  LPj.  The  basis  and  dual  variables  are  also 
optimal  for  LPo,  and  the  required  solution  is  x  =  xq.  The  procedure  may 
be  terminated  after  x  and  its  bounds  are  restored  to  their  original  values  (xo 
and  /). 

2.  LPi  is  found  to  be  unbounded  when  a  certain  nonbasic  variable  is  considered 
for  entry  into  the  basis.  The  same  basis  and  nonbasic  variable  will  produce  a 
feasible  descent  direction  for  LPo-  (This  is  the  direction  of  recession  described 
by  Osborne  [Osb85].)  The  variables  and  bounds  are  again  restored  to  their 
original  values,  and  iterations  continue  on  the  original  problem. 

3.  A  degenerate  vertex  arises  (say  x  =  xj).  We  may  again  invoke  the  procedure 
of  defining  a  subsidiary  LP.  Since  the  variable  that  just  entered  the  basis  must 
have  moved  away  from  its  bound  in  order  to  cause  the  degeneracy,  the  degree 
of  degeneracy  must  be  less  than  before.  The  procedure  can  be  invoked  again, 
using  x\  to  define  a  new  subsidiary  problem  LP2. 

In  general,  the  procedure  may  be  applied  recursively.  Starting  with  k  —  0, 
problem  LP*  reaches  case  1,  2  or  3  in  a  finite  number  of  iterations  (since  the  objective 
function  for  LP*  is  monotonically  decreasing).  Case  3  leads  to  a  new  problem 
LPfc+i  but  can  occur  only  a  finite  number  of  times  (since  the  degree  of  degeneracy 
is  monotonically  decreasing). 

0.2.  Discussion 

Wolfe’s  procedure  is  appealing  for  at  least  two  reasons:  it  uses  the  simplex  method 
itself  to  resolve  degeneracy,  and  it  can  be  implemented  with  a  minimum  of  overhead 
(at  least  for  the  case  l  =  0,  u  =  00),  as  shown  by  Ryan  and  Osborne  [RO86]. 

Although  a  degenerate  vertex  is  unlikely  to  be  encountered  in  LP*  ( k  >  0), 
particularly  if  d  is  defined  using  random  positive  numbers,  it  remains  necessary  to 
cope  with  the  possibility.  This  may  be  an  inconvenience. 

Another  drawback  is  the  need  to  decide  that  degeneracy  is  present  and  the 
need  to  define  the  precise  set  of  degenerate  constraints.  Suppose  we  have  a  set  of 
basic  variables  that  are  not  quite  on  their  bounds.  This  could  include  all  the  basic 
variables.  If  we  do  not  include  some  of  them  in  the  definition  of  xD,  it  is  probable 
that  only  a  very  short  step  will  be  taken  after  the  degeneracy  has  been  “resolved”, 
and  we  may  need  to  invoke  the  procedure  again. 

6.3.  A  modification 

Suppose  we  introduce  a  parameter  6  into  the  Wolfe  procedure,  so  that  instead  of 
xD  >  l  —  d  we  now  have  xD  >  l  —  6d,  where  6  >  0.  We  shall  regard  d  as  fixed,  but  5 
remains  to  be  specified. 

Note  that  if  there  is  a  unique  choice  of  blocking  variable  when  6  =  I,  the  same 
blocking  variable  will  be  chosen  for  any  positive  value  of  6.  Even  if  the  choice  is 
not  unique,  providing  we  use  a  consistent  criterion  for  choosing  among  the  set  of 
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blocking  variables,  the  choice  will  still  be  independent  of  6.  (Either  of  the  EXPAND 
procedures  would  be  suitable  for  making  the  choice.) 

We  may  extend  the  argument  to  show  that  the  sequence  of  bases  generated  when 
solving  LPi  is  independent  of  6. 

The  significance  of  this  observation  is  that  6  may  be  chosen  extremely  small 
(assuming  ||dj|  is  of  order  1).  But  then,  if  6  is  sufficiently  small  there  is  no  need  to 
define  LPj;  we  can  simply  solve  the  original  problem  with  a  tiny  modification  to 
some  of  the  bounds. 

We  emphasize  that  solving  LPo  with  the  modified  bounds  generates  the  same 
bases  as  solving  LPi,  and  if  no  further  degenerate  vertices  are  encountered,  we  either 
determine  a  feasible  descent  direction  or  confirm  that  x  is  the  required  solution. 


6.4.  A  further  modification 

As  in  the  original  Wolfe  procedure,  the  difficulty  is  that  we  cannot  guarantee  that 
a  further  degenerate  vertex  will  not  occur.  We  now  show  how  to  avoid  such  an 
eventuality  by  judicious  choice  of  d.  First  note,  however,  that  choosing  the  elements 
of  d  at  random  is  not  a  good  strategy  from  the  perspective  of  preserving  well- 
conditioned  bases.  Indeed,  the  best  choice  is  d  =  e,  the  vector  of  ones  (assuming 
the  problem  is  well  scaled).  In  this  case,  the  blocking  variable  corresponds  to  the 
largest  eligible  pivot  element  in  the  search  direction  p. 

Such  a  structured  choice  for  d  would  seem  to  increase  the  probability  of  a  degen¬ 
erate  vertex  arising.  Observe,  however,  that  once  the  best  possible  choice  of  blocking 
variable  has  been  defined  (say  xr),  only  dT  need  be  defined.  We  are  then  free  to 
increase  d,,  i  ^  r,  since  had  this  been  our  initial  choice  of  d,  the  blocking  variable 
would  remain  the  same.  Suppose  we  increase  all  elements  d;  (t  /  r)  by  r,  where  r 
is  of  order  one.  Such  a  choice  prevents  the  current  vertex  from  being  degenerate  in 
LPj.  Also,  in  the  following  iteration  the  value  of  the  step  to  the  blocking  variable  is 
relatively  large  (of  order  r),  and  hence  again  results  in  a  good  choice  for  the  blocking 
variable  as  regards  the  condition  of  the  basis. 

The  second  iteration  fixes  dr<,  say,  but  we  are  still  free  to  alter  d,,  t  ^  r,  r' . 
By  proceeding  in  this  manner  we  can  prevent  the  occurrence  of  degenerate  vertices. 
Note  that  although  dr  has  been  fixed,  if  xr  ever  becomes  basic  again,  we  are  able  to 
redefine  its  bound.  This  may  appear  to  negate  the  argument  that  it  would  still  be 
the  first  blocking  variable.  However,  provided  all  di  are  being  increased  similarly, 
there  is  no  contradiction. 

With  6  small  and  d  =  6e  used  to  modify  all  bounds  (not  just  those  of  xD),  the 
approach  just  described  is  equivalent  to  the  EXPAND  procedure. 

Note  that  a  strict  implementation  of  Wolfe’s  procedure  would  preserve  x0  and 
eventually  determine  a  feasible  descent  direction  p.  In  the  modified  procedure,  we 
take  a  sequence  of  steps  (to  x  say)  before  p  is  determined.  We  then  step  along  p 
from  x  rather  than  xo-  However,  if  6  is  sufficiently  small,  ||x0  -  x|[  is  negligible  and 
the  infeasibility  incurred  (with  respect  to  x  >  l)  is  strictly  bounded. 
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6.5.  Summary 

We  have  shown  that  the  EXPAND  procedure  is  closely  related  to  Wolfe’s  method. 

Some  advantages  are  as  follows: 

•  There  is  no  need  to  judge  whether  or  not  degeneracy  is  present,  or  to  specify 
a  set  of  degenerate  variables. 

•  There  is  no  storage  overhead  or  logical  overhead.  General  bounds  on  x  can  be 
handled  without  complication. 

•  There  is  no  numerical  or  logical  information  to  be  preserved  across  basis  fac¬ 
torizations  (other  than  x  and  the  current  6). 

•  All  iterations  are  “equal”,  in  the  sense  that  there  is  exactly  one  steplength 
determination  per  iteration.  (In  Wolfe’s  method,  if  LP*  is  found  to  be  degen¬ 
erate,  the  steplength  procedure  is  effectively  repeated  at  the  beginning  and 
end  of  LPfc+i.) 


7.  Issues  Arising  in  Phase  1 

Broadly  speaking,  Phase  1  of  the  simplex  method  is  implemented  by  applying  a 
normal  (Phase  2)  procedure  to  a  modified  problem ,  whose  bounds  have  been  altered 
to  make  the  current  solution  feasible. 

First  we  need  to  summarize  the  main  aspects  of  Phase  1.  (See  also  Orchard-Hays 
[Orc68]  and  Beale  [Bea70].)  Some  finer  points  can  then  be  discussed. 

7.1.  The  Phase- 1  bounds  and  objective 

Suppose  a  variable  Xj  has  true  bounds  ( lj,uj ).  If  xj  lies  above  its  upper  bound  by 
more  than  the  current  feasibility  tolerance  ( xj  —  Uj  >  6k),  its  bounds  are  treated  as 
(Zj,oo).  Similarly,  if  its  lower  bound  is  not  satisfied  ( lj  -  Xj  >  6k),  the  bounds  are 
taken  to  be  (-00 ,Uj).  Otherwise,  the  true  bounds  on  Xj  are  retained. 

The  Phase- 1  objective  function  cTx  is  then  defined  as  Cj  =  1,  -1  and  0  re¬ 
spectively.  Note  that  c  is  redefined  every  Phase-1  iteration.  It  is  used  to  compute 
reduced  costs  and  hence  a  search  direction  p  in  the  normal  way. 

For  most  Phase- 1  iterations,  an  appropriate  step  along  p  can  be  defined  as  usual 
by  applying  the  EXPAND  procedure  to  the  data  ( x,pj,u ),  where  Z  and  u  are  the 
modified  bounds.  We  define  this  step  to  be  op,  since  it  keeps  variables  feasible  with 
respect  to  any  bounds  that  they  currently  satisfy. 

For  some  iterations,  a  special  step  a/  must  be  taken;  see  Section  7.3.  (This  step 
allows  one  or  more  infeasible  variables  to  become  feasible.) 

Thus,  Phase  1  is  essentially  the  same  as  Phase  2  except  that  c,  Z  and  u  are 
redefined  every  iteration,  and  two  possible  steps  are  computed  rather  than  one. 
(The  step  ultimately  taken  is  a  =  ap  or  a/,  subject  to  safeguards  discussed  below.) 

To  see  that  progress  is  guaranteed,  observe  that  the  Phase- 1  objective  decreases 
as  a  increases  from  zero.  If  any  infeasible  variables  become  feasible  as  a  continues 
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to  increase,  the  sum  of  infeasibilities  will  decrease  at  a  lower  rate  (and  could  even 
start  to  increase),  but  the  number  of  infeasibilities  will  be  lower.  Convergence  is 
therefore  assured. 

(A  more  intricate  steplength  procedure  can  be  designed  to  minimize  the  piece- 
wise  linear  function  cT(x  +  ap),  where  c  is  regarded  as  a  function  of  a;  for  example, 
see  Greenberg  [Gre78],  Fourer  [Fou85j.  However,  we  adopt  the  simpler  approach  as 
it  is  effective  in  practice.  Both  approaches  have  the  desirable  property  that  many 
infeasibilities  can  be  removed  in  one  iteration.) 

7.2.  Advantages  of  increasing  8k 

Note  that  the  EXPAND  procedure  does  allow  a  to  increase  from  zero.  Also,  if  the 
final  steplength  leaves  the  number  of  infeasibilities  unaltered,  the  sum  of  infeasi¬ 
bilities  at  the  start  of  the  next  Phase- 1  iteration  will  be  lower  simply  because  the 
tolerance  used  to  measure  infeasibility  has  increased  (from  8k  to  6k  +  r).  Thus  for 
two  separate  reasons,  either  the  sum  or  the  number  of  infeasibilities  will  decrease 
after  each  Phase- 1  iteration. 

Although  r  is  typically  very  small,  it  is  intended  to  be  significantly  larger  than 
machine  precision  e,  and  preferably  larger  than  the  cut-olT  value  toip  (Section  3.1). 
An  important  benefit  is  that  it  helps  mask  the  rounding  error  that  is  inevitably 
present  when  x  is  updated  to  x  +  ap.  The  set  of  infeasible  variables  (as  measured  by 
6k)  is  therefore  guaranteed  to  stay  the  same  or  to  diminish.6  We  believe  that  many 
“infinite  loop”  failures  of  simplex  implementations  have  been  due  to  an  inadvertent 
oscillation  in  the  number  of  infeasibilities  when  c  is  redefined  each  Phase- 1  iteration 
with  a  constant  6k-  (An  example  is  described  by  Ogryczak  [Ogr87].  Similar  examples 
were  encountered  with  MINOS  prior  to  the  present  implementation.) 

If  a  problem  is  still  infeasible  after  K  iterations,  the  feasibility  tolerance  is  re¬ 
duced  from  6k  to  60  for  the  next  cycle  of  iterations.  An  apparent  disadvantage  is  that 
the  number  of  infeasibilities  may  increase  by  some  arbitrary  number  (say  q),  and  the 
sum  could  increase  by  as  much  as  qbK-  However,  this  is  normally  inconsequential 
even  if  q  is  nearly  as  large  as  m.  (Here  it  is  important  that  many  infeasibilities  can 
be  removed  in  one  Phase-1  iteration.)  Similar  comments  apply  when  the  resetting 
procedure  is  invoked  at  an  apparently  optimal  solution. 

7.3.  The  special  Phase- 1  step 

By  construction,  the  Phase- 1  objective  causes  at  least  some  of  the  infeasible  variables 
to  move  towards  the  feasible  region.  Let  dj  denote  the  step  that  allows  such  a 
variable  to  reach  its  nearest  bound  exactly  (and  hence  become  feasible): 


(  (h  -  *j)/pj 


Pi  >  0, 

Pi  <  0, 
otherwise 


8Fletcher’s  method  for  resolving  degeneracy  also  possesses  favorable  prop 
of  rounding  error;  see  [Fle85,Fle87], 
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(cf.  (6)).  The  special  Phase- 1  step  referred  to  above  is  then  a/  —  d:  for  some  index 
j  G  I  (where  /  temporarily  denotes  the  set  of  infeasible  variables). 

In  general  it  would  be  sensible  to  define 


a/ =  a,  =  majcdj,  a  =  min{or/,ap},  (12) 

since  if  any  infeasible  variables  become  feasible  as  the  steplength  increases,  aj  marks 
the  point  at  which  the  maximum  number  become  feasible.  (Note  that  a  should  not 
exceed  a/  because  some  of  the  infeasible  variables  could  become  more  infeasible  as 
a  increases.  Also  note  that  ap  could  be  infinite.) 

If  a/  =  a,,  a  danger  is  that  the  corresponding  pivot  element  p,  could  be  arbi¬ 
trarily  small. 

Following  the  philosophy  of  Harris  (Section  3.2),  some  freedom  to  maximize  the 
pivot  element  can  be  obtained  by  using  a  two-pass  procedure.  The  perturbed  bounds 
l  -  6e  and  u  +  6e  are  used  in  the  first  pass  as  usual,  but  this  gives  a  step  al  that 
is  slightly  too  small.  The  second  pass  then  considers  all  unperturbed  steps  <5,  no 
smaller  than  al: 

|p,|  =  max|pj|  such  that  dj>al. 

This  was  the  method  used  in  Lustig’s  experiments  [Lus87]  and  in  all  preceding 
versions  of  MINOS.  It  appears  to  have  performed  reliably  for  many  years. 

Nevertheless  it  is  evident  that  |p,|  could  still  be  as  small  as  the  cut-off  value  tolp. 
We  have  therefore  adopted  the  following  safer  stategy.  In  the  first  pass,  we  find  the 
largest  relevant  pivot  element: 

<t>  =  (13) 

In  the  second  pass  we  then  find  the  largest  step  subject  to  the  pivot  element  being 
reasonably  close  to  <f>: 

ar  =  a,  =  maxdj  such  that  |p>(  >  -f<t>  (14) 

for  some  constant  7,  where  0  <  7  <  1.  Experience  suggests  that  the  step  a /  should 
be  taken  whenever  possible  (to  remove  the  associated  infeasible  variable  from  the 
basis).  We  therefore  define 


»-{£  Zi:L, 

where  al  (from  the  first  pass  of  the  Harris  procedure)  is  slightly  larger  than  ap. 
Together,  ( 13)— ( 15)  define  the  steplength  in  place  of  (12).  By  observation  on  53  test 
problems,  the  values  7  =  0.1  and  7  =  0.01  seem  to  impede  Phase  1  relative  to  the 
unsafeguarded  7  =  0.  We  have  therefore  settled  on  7  =  0.001. 

A  final  comment:  computation  of  the  special  step  a  j  does  not  depend  criti¬ 
cally  on  the  feasibility  tolerance,  and  is  therefore  compatible  with  the  EXPAND 
procedure. 


so 


An  anti-cycling  procedure 


8.  Nonlinear  Programs  with  Linear  Constraints 

Here  we  assume  that  the  problem  to  be  solved  is  the  same  as  in  (1),  except  that  the 
objective  function  cTx  is  replaced  by  a  smooth  general  function  F(x).  The  problem 
could  therefore  be  a  quadratic  program  (QP)  or  a  more  general  linearly  constrained 
(LC)  optimization  problem. 

The  environment  we  have  in  mind  is  active-set  methods  for  QP  and  LC  prob¬ 
lems  (Fletcher  [Fle81];  Gill,  Murray  and  Wright  [GMW81]).7  It  has  been  observed 
by  Osborne  [Osb85]  that  Wolfe’s  anti-cycling  procedure  generalizes  to  certain  LC 
algorithms,  including  the  reduced-gradient  method  of  Wolfe  [Wol62].  The  same  is 
true  of  the  present  approach.  The  following  preliminary  implementation  has  been 
developed  for  the  reduced-gradient  algorithm  in  MINOS. 

8.1.  A  normal  iteration 

Conceptually,  the  EXPAND  procedure  may  be  applied  directly.  Thus,  at  each  it¬ 
eration  we  increase  the  feasibility  tolerance  slightly,  and  after  obtaining  a  search 
direction  p,  we  compute  a  positive  step  and  a  blocking  index  (a,  r)  as  before.  Wc 
rename  this  step  qm«»  since  it  may  be  preferable  to  take  a  shorter  step.  Any  step 
in  the  interval  f0.qm«,]  will  give  a  point  that  is  acceptably  feasible. 

In  general,  we  then  perform  a  linesearch  to  find  a  step  cr  that  approximately 
minimizes  the  objective  function  over  the  specified  interval: 

d  ss  argmin  F{x  +  op),  a  G  (O,^]. 


Note  that 

Qm.v  — -  max{oaun,q2}, 

where  a 2  is  the  step  that  allows  the  blocking  variable  to  reach  its  true  bound  (see 
Figure  2).  In  case  1  of  Figure  2,  the  search  would  be  performed  over  the  “large” 
interval  (0,o2),  while  for  cases  2  and  3  it  would  be  performed  over  (0,q,nin]. 

If  the  linesearch  is  successful  (the  objective  function  is  “sufficiently  reduced”), 
there  is  no  danger  of  cycling  and  the  optimization  algorithm  proceeds  normally.  The 
current  point  is  updated  (x  <—  i  •+■  dp),  and  if  the  maximum  step  was  chosen,  the 
blocking  constraint  is  added  to  the  working  set  (i.e.,  xr  becomes  nonbasic). 

8.2.  Avoiding  the  linesearch 

In  practice,  it  may  be  inefficient  or  unwise  to  attempt  a  linesearch,  since  Omax  could 
be  very  small  if  degeneracy  is  present.  Even  if  the  linesearch  returns  the  maximum 
step  d  =  Omax,  the  improvement  in  objective  value  may  be  very  slight.  More 
seriously,  the  “noise  level”  in  F(x  +  ap)  over  the  interval  (0.  o™.»]  may  be  too  great 
to  allow  identification  of  an  improved  point,  and  the  linesearch  will  be  obliged  to 
“fail”. 


rThe  EXPAND  procedure  ha*  been  implemented  in  the  1988  versions  of  QPSOL  and  LSSOL. 
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We  therefore  make  use  of  the  step  a2  that  allows  the  blocking  variable  to  reach 
its  bound  exactly.  In  Figure  2,  this  step  is  positive  for  the  first  two  cases,  but 
negative  for  the  third  (and  therefore  not  shown). 

In  cases  1  and  2  ( a2  >  0)  we  always  perform  the  linesearch. 

In  case  3  ( a2  <  0),  degeneracy  is  present  and  we  usually  take  a  zero  step.  (We 
skip  the  linesearch  and  make  the  blocking  variable  nonbasic.) 

The  only  exception  is  when  case  3  would  create  a  vertex  of  the  feasible  region; 
this  is  the  most  likely  circumstance  under  which  a  zero  step  could  lead  to  cycling. 
If  a2  <  0  and  the  working  set  has  only  one  degree  of  freedom,  we  attempt  to  find  a 
positive  step  by  performing  a  linesearch  over  the  interval  (0,  Qmin]- 

8.3.  Recovery  from  a  linesearch  failure 

If  the  linesearch  ever  fails  to  find  an  improved  point,  we  try  to  determine  whether 
the  failure  was  due  to  the  search  interval  being  too  small. 

If  a2  <  Omin  (cases  2  and  3),  we  force  a  step  to  the  true  bound  (a  =  a2)  and 
update  the  working  set. 

Otherwise,  we  assume  that  a  better  search  direction  is  required.  We  leave  the 
working  set  unaltered  and  invoke  a  series  of  recovery  procedures.  (In  MINOS,  these 
include  re-estimating  any  unknown  components  of  the  objective  gradient  using  cen¬ 
tral  differences,  resetting  the  approximate  reduced  Hessian,  deleting  a  constraint 
from  the  working  set,  and  refactorizing  the  basis.) 

8.4.  Nonlinear  constraints 

Nonlinear  programs  involving  nonlinear  constraints  are  often  treated  by  SQP  and 
SLC  methods,  involving  a  sequence  of  linearly  constrained  subproblems  to  which 
the  above  anti-degeneracy  procedures  may  be  applied. 

It  is  clearly  necessary  that  cycling  be  avoided  within  the  LC  subproblems. 
Whether  this  is  sufficient  to  prevent  the  overall  NLC  algorithm  from  cycling  re¬ 
mains  to  be  investigated. 
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9.  Computational  Results 

In  this  section  we  compare  three  steplength  procedures  for  the  simplex  method.  For 
convenience  we  give  them  the  following  names: 


SP1:  The  “textbook”  ratio  test  of  Section  3.1. 

SP2:  The  simplified  EXPAND  procedure  of  Section  5. 

SP3:  The  maximum-pivot  EXPAND  procedure  of  Section  4.  (This  includes  Harris- 
type  tie-breaking  and  is  the  preferred  method  of  this  paper.) 

SP3  has  been  implemented  in  GAMS/MINOS  (see  [BKM88])  and  in  MINOS  5.3 
(May  1988).  All  three  procedures  have  been  implemented  in  MINOS  5.3  and  tested 
on  the  first  53  linear  programs  in  the  netlib  collection  [Gay 85].  The  problems  were 
ordered  according  to  the  number  of  nonzero  elements  as  in  [Lus87].  The  main  run- 


time  options  specified  were 

PRINT  LEVEL 

0 

CRASH  0PTI0I 

1 

CRASH  T0LERAICE 

0.1 

SCALE  OPTION 

2 

PARTIAL  PRICE 

10 

EXPAND  FREQUENCY 

10000 

FEASIBILITY  TOLERANCE 

1.0E-6 

(the  default  options  for  linear  problems  in  MINOS  5.3).  The  last  two  options  define 
I(  =  10000  and  8  =  10~6  for  the  EXPAND  procedures.  The  limit  on  calls  to  the 
resetting  procedure  at  an  apparent  optimum  was  set  to  ft  =  2  (Section  4.4). 

The  CRASH  parameters  cause  MINOS  to  choose  an  approximately  triangular  basis 
from  the  columns  of  A  =  (  A  I).  In  most  cases  the  scaling  option  has  the  effect  of 
making  ||x*||  =  0(1),  where  x*  is  the  scaled  optimal  solution.8  This  helps  justify 
8  =  10-6  as  a  feasibility  tolerance  for  the  scaled  problem. 

Tables  1,  2  and  3  give  results  for  the  simplex  method  using  SP1,  SP2  and  SP3 
respectively.  The  “objective  function”  values  indicate  that  the  final  objective  was 
accurate  to  four  or  more  digits  (except  for  one  problem  that  terminated  early).  The 
meaning  of  “degenerate  steps”  depends  on  the  method;  see  below.  Solution  times  are 
given  in  CPU  seconds;  they  do  not  include  time  for  data  input  or  solution  output.9 

Figure  4  plots  the  “total  iterations”  in  Tables  1  and  2  relative  to  the  iterations 
in  Table  3.  Figure  5  compares  CPU  times  similarly. 


’Exceptions  were  problems  GROW7,  GROW15  and  GROW22,  for  which  ||x*||  =  O(107),  ||x*||  = 
0(10*). 

’Tests  were  run  as  batch  jobs  on  a  DEC  VAXstation  II.  The  operating  system  was  VAX/VMS 
version  4.5.  The  compiler  was  VAX  FORTRAN  version  4.6  with  default  options,  including  code 
optimization  and  DJloating  arithmetic  (relative  precision  r  ss  2.8  x  10-17).  The  memory  available 
kept  paging  to  a  minimum. 


9.  Computational  Results 
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9.1.  The  textbook  ratio  test 

SP1  was  safeguarded  by  treating  pj  as  zero  in  equation  (6),  using  tolp  =  c2/3.  Since 
rounding  error  can  cause  the  steplength  to  be  negative,  a  further  precaution  was  to 
set  a  =  0  if  the  ratio  test  gave  a  <  10-16.  (“Degenerate  steps”  counts  the  number 
of  times  this  occurred.)  Following  conventional  practice,  blocking  variables  were  set 
exactly  on  their  bounds  when  they  became  nonbasic. 

Although  it  would  be  reasonably  easy  to  break  (near)  ties  in  favor  of  large  | , 
we  chose  not  to  tamper  further  with  the  classical  procedure;  methodical  tie-breaking 
is  the  province  of  the  Harris  and  EXPAND  procedures. 

In  the  test  runs,  small  pivots  slipped  through  the  c2/3  sieve  several  times  on  each 
of  five  problems  (SCTAP2,  SCTAP3,  SCSD8,  PILOTJA  and  PILOT).  In  general,  these 
are  detected  as  near-singularities  when  the  LU  factors  of  the  basis  are  updated. 
Refactorization  is  invoked  and  some  variable  Xj  is  replaced  by  an  appropriate  slack 
variable.  Since  x}  retains  its  value  when  rejected  from  the  basis,  iterations  continue 
without  apparent  interruption. 

Only  one  real  failure  was  encountered:  problem  SCSD8  was  terminated  at  itera¬ 
tion  5804  after  stalling  for  the  final  1000  iterations  (far  short  of  optimality).  There 
was  no  obvious  cycle  in  progress,  but  small  pivots  were  encountered  frequently  dur¬ 
ing  the  run,  causing  the  basis  to  be  ill-conditioned  for  many  groups  of  iterations. 
Empirically,  ill-conditioning  can  only  aggravate  stalling  (particularly  for  a  method 
that  has  no  guarantee  of  terminating). 

Figure  4  illustrates  that  on  most  problems,  SP1  led  to  more  simplex  iterations 
than  SP2  or  SP3. 

9.2.  The  simplified  EXPAND  procedure 

In  Table  2,  “degenerate  steps”  means  the  number  of  times  ST2  required  two  passes 
to  determine  a  blocking  variable. 

No  singularities  were  encountered  during  the  tests,  and  all  problems  terminated 
successfully.  The  resetting  procedure  was  invoked  at  10000  iterations  for  the  last 
two  problems,  with  372  and  423  nonbasic  variables  (respectively)  being  moved  a 
distance  in  the  range  (c08,^-)  onto  their  bounds. 

With  6  as  small  as  10-6,  resets  do  not  disturb  x  greatly.  After  resetting  at  an  ap¬ 
parent  optimum,  most  problems  were  confirmed  optimal  with  no  further  iterations. 
On  PILOT4,  GANGES,  PILOTJA  and  PILOT,  35,  190,  98  and  28  nonbasic  variables 
were  moved  onto  their  bound  and  4,  1,  1  and  15  additional  iterations  were  needed 
to  confirm  optimality.  In  the  case  of  GANGES,  a  second  reset  moved  one  nonbasic 
variable  but  led  to  no  further  iterations. 


9.3.  The  EXPAND  procedure 

In  Table  3,  “degenerate  steps”  means  the  number  of  times  a  was  forced  to  be  as 
large  as  amjn.  This  is  the  number  of  times  a  blocking  variable  was  made  nonbasic 
at  an  infeasible  value,  rather  than  reaching  its  bound  exactly. 
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An  anti-cycling  procedure 


To  date,  there  has  been  no  specific  study  of  the  effect  of  maximizing  the  pivot 
e’ement  within  a  steplength  procedure — the  Harris  approach  to  tie-breaking.  Folk¬ 
lore  has  it  that  “stability  is  improved  and  the  number  of  simplex  iterations  is  often 
reduced”.  However,  such  a  statement  is  not  especially  meaningful  without  a  precise 
definition  of  the  procedures  being  compared. 

Here,  it  is  meaningful  to  compare  both  versions  of  the  EXPAND  procedure  be¬ 
cause  it  is  understood  that  in  each  case  the  surrounding  simplex  algorithm  increases 
the  feasibility  tolerance  every  iteration  and  deals  correctly  with  infeasible  blocking 
variables  (by  retaining  their  numerical  values  when  they  become  nonbasic). 

The  results  in  Tables  2  and  3  essentially  confirm  the  folklore.  Figure  4  illustrates 
the  trend  more  clearly:  the  ratio  of  SP2  to  SP3  iterations  is  mostly  greater  than 
one.  In  Figure  5,  the  CPU-time  ratios  are  shifted  slightly  downwards,  reflecting  the 
fact  that  SP2  usually  requires  only  one  pass,  whereas  SP3  always  requires  two. 

Only  two  problems  required  additional  iterations  after  resetting.  On  PILOTJA, 
81  nonbasics  were  movrd  onto  their  bound  at  iteration  6070,  and  63  iterations  later 
a  second  reset  moved  7  nonbasics.  Termination  occurred  8  iterations  later  (with  no 
attempt  to  reset). 

On  PILOT,  127  nonbasics  were  moved  onto  their  bounds  at  iteration  10000, 
and  43  at  iteration  17698.  Termination  occurred  after  a  further  18  iterations.  As 
expected,  the  effect  of  resets  was  slightly  less  noticeable  than  with  SP2. 

9.4.  Other  parameter  values 

The  53  test  problems  have  been  solved  many  times,  with  and  without  scaling  and 
partial  pricing.  One  of  the  main  parameters  of  interest  is  the  feasibility  tolerance. 
We  have  experimented  with  the  values  6  -■  1 0— 4 ,  10-5,  10-6  and  10~7  (Harris 
recommended  6  =  5  X  10~4),  but  the  sensitivity  of  the  simplex  method  to  minor 
algorithmic  changes  seems  to  have  masked  any  useful  trend.  Significant  improve¬ 
ments  were  certainly  observed  on  some  of  the  problems  with  6  =  10~4.  The  risk  is 
a  greater  disturbance  aft^r  resetting  on  problems  that  are  somewhat  ill-conditioned 
(notably  PILOTJA  and  PILOT). 

As  a  further  test,  we  disabled  the  “expand”  feature  of  SP3  by  specifying  A'  =  oo, 
t  =  0.  This  has  the  effect  of  fixing  the  feasibility  tolerance  at  60  =  0.5  x  10~6,  and 
most  closely  resembles  the  Harris-type  tie-breaking  (without  losing  the  feature  of 
making  infeasible  blocking  variables  nonbasic  at  their  correct  value).  No  failures 
occurred  on  four  runs  with  and  without  scaling  and  partial  pricing,  confirming  that 
the  probability  of  failure  with  the  Harris  procedure  is  indeed  low,  given  the  second 
safeguard.  The  iteration  counts  were  much  the  same  as  when  6  was  allowed  to 
expand. 

Note  that  once  the  second  safeguard  is  implemented,  the  assurance  gained  by 
allowing  6  to  expand  comes  at  no  cost.  There  is  no  reason  to  keep  6  fixed. 


9.  Computational  Results 
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Problem  Objective  Total  Degen  Percent  Solve  tine 

function  itns  steps  degen  VAX  II  secs 


1 

AFIR0 

-4 . 64753 14285714E+02 

6 

3 

50.00 

0.45 

2 

ADLITTLE 

2. 25494963 16238E+05 

113 

20 

.7.70 

5.90 

3 

SC205 

-5 . 220206121 1707E+01 

118 

14 

11.86 

13.03 

4 

SCAGR7 

-2. 3313897523795E+06 

86 

10 

11.63 

6.80 

S 

SHARE 2B 

-4. 1673224074142E+02 

124 

48 

38.71 

8.31 

6 

RECIPE 

-2 . 6661600000000E+02 

33 

3 

9.09 

2.12 

7 

VTPBASE 

1.2983146246136E+05 

64 

30 

46.88 

6.39 

8 

SHARE1B 

-7. 658931 85791 86E +04 

277 

4 

1.44 

24.68 

9 

B0RE3D 

1 . 3730803942085E+03 

207 

155 

74.88 

28.25 

10 

SCORPIOS 

1.8781248227381E+03 

109 

41 

37.61 

21.23 

11 

CAPRI 

2 . 6900 1 29 1 37682E+03 

246 

50 

20.33 

32.86 

12 

SCAGR25 

-1 . 4753433060769E+07 

376 

79 

21.01 

89.30 

13 

SCTAP1 

1.4122500000000E+03 

259 

123 

47.49 

38.68 

14 

BRANDY 

1 . 5185098964881E+03 

331 

74 

22.36 

53.65 

15 

ISRAEL 

-8 . 9664482 186305E+05 

266 

20 

7.52 

35.55 

16 

ETAMACR0 

-7 . 557 152 1755 166E+02 

470 

134 

28.51 

93.29 

17 

SCFXM1 

1 . 8416759028349E+04 

378 

86 

22.75 

69.71 

18 

GR0V7 

-4.7787811814712E+07 

179 

65 

36.31 

37.90 

19 

BANDM 

- 1.5862801 8450 12E+02 

483 

80 

16.56 

100.91 

20 

E226 

-1 . 8751929066371E+01 

580 

202 

34.83 

85.34 

21 

STANDATA 

1 . 2576995000000E+03 

145 

109 

75.17 

23.71 

22 

SCSD1 

8.6666666743334E+00 

1084 

1000 

92.25 

89.14 

23 

GFRDPNC 

6 . 9022359995488E+06 

672 

350 

52.08 

187.48 

24 

BEAC0NFD 

3 . 3592485807200E+04 

116 

27 

23.28 

13.90 

25 

STAIR 

-2. 5126695119296E+02 

469 

63 

13.43 

231.65 

26 

SCRS8 

9. 0429998618888E+02 

537 

183 

34.08 

143.88 

27 

SEBA 

1.571 1600000000E+04 

445 

54 

12.13 

106.98 

28 

SHELL 

1 . 2088253460000E+09 

303 

73 

24.09 

76.90 

29 

PIL0T4 

-2.581 1392588836E+03 

1613 

229 

14.20 

706.74 

30 

SCFXH2 

3 . 6660261564999E+04 

917 

201 

21.92 

326.56 

31 

SCSD6 

5 . 0500000078262E+01 

1476 

1137 

77.03 

214.56 

32 

GROW 15 

- 1.0687094 129358E+08 

397 

148 

37.28 

150.31 

33 

SHIP04S 

1 . 7987 147  004454E+06 

152 

34 

22.37 

36.23 

34 

FFFFF800 

6 . 5567959102690E+05 

938 

388 

41.36 

287.07 

35 

GANGES 

-1 . 0958596920679E+05 

687 

214 

31.15 

358.13 

36 

SCFXH3 

5. 490 125454975 1E+04 

1359 

303 

22.30 

724.11 

37 

SCTAP2 

1 . 7248071428571E+03 

785 

630 

67.52 

453.01 

38 

GR0U22 

- 1 . 6083433648256E+08 

635 

254 

40.00 

371.58 

39 

SHIP04L 

1.7933245379704E+06 

266 

55 

20.68 

65.46 

40 

PILOTWE 

-2 . 7200970057530E+06 

5003 

838 

16.75 

3540.04 

41 

SIERRA 

1 . 5394460531792E+07 

1340 

825 

61.57 

704.48 

42 

SHIP08S 

1 . 9200982105346E+06 

242 

64 

26.45 

102.41 

43 

SCTAP3 

1 . 42400000 00000E +03 

1151 

876 

76.11 

837.87 

44 

SHIP 12S 

1.4892361344061E+06 

434 

111 

25.58 

253.79 

45 

25FV47 

5. 5018458882865E+03 

8687 

1100 

12.66 

5794.72 

46 

SCSD8 

1 . 318E+03  Stalled 

5804 

3566 

61.44 

2787.51 

47 

NESM 

1. 4076079386 175E+07 

3067 

0 

0.00 

1338.22 

48 

CZPROB 

2. 1851966988566E+06 

1519 

151 

9.94 

788.61 

49 

PILOT JA 

-6. 1130625520046E+03 

7086 

504 

7.11 

6317.42 

50 

SHIP08L 

1. 90905521 13891E+06 

522 

149 

28.54 

250.04 

51 

SHIP12L 

1.4701879193293E+06 

877 

252 

28.73 

594.19 

52 

80BAU3B 

9 . 87230019 10483E+05 

10896 

1483 

13.61 

12674.42 

53 

PILOT 

-5 . 5740380062649E+02 

20720 

1630 

7.87 

81173.13 

Table  1:  Results  from  MINOS  5.3  using  textbook  ratio  test 
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An  anti-cycling  procedure. 


Problem 

Objective 

function 

Total 

itns 

Degen 

steps 

Percent 

degen 

Solve  time 

VAX  II  secs 

1 

AFIRO 

-4 . 64753 142857 14E+02 

6 

3 

50.00 

0.47 

2 

ADLITTLE 

2. 2549496316238E+05 

92 

12 

13.04 

4.75 

3 

SC205 

-5 . 220206 121 1707E+0 1 

123 

13 

10.57 

14.23 

4 

SCAGR7 

-2. 3313897523795E+06 

86 

10 

11.63 

6.73 

5 

SHARE2B 

-4. 1S73224074142E+02 

113 

17 

15.04 

7.69 

6 

RECIPE 

-2.666 1600000000E+02 

33 

3 

9.09 

1.91 

7 

VTPBASE 

1.29831 46246 136E+05 

71 

24 

33.80 

6.97 

8 

SHARE IB 

-7.658931 857 9 1 86E+04 

242 

3 

1.24 

21.56 

9 

B0RE3D 

1.3730803942085E+03 

165 

80 

48.48 

23.41 

10 

SC0RPI0H 

1 . 8781248227381E+03 

105 

41 

39.05 

20.10 

11 

CAPRI 

2 . 6900129137682E+03 

245 

22 

8.98 

33.21 

12 

SCAGR25 

-1 . 4753433060769E+07 

361 

68 

18.84 

87.47 

13 

SCTAP1 

1.4122500000000E+03 

291 

70 

24.05 

42.55 

14 

BRANDY 

1 . 5185098964881E+03 

462 

35 

7.58 

75.15 

15 

ISRAEL 

-8 . 9664482 186305E+05 

225 

29 

12.89 

30.13 

16 

ETAHACRO 

-7 . 557 15217 18413E+02 

570 

135 

23.68 

119.26 

17 

SCFXM1 

1 . 8416759028349E+04 

396 

61 

15.40 

72.28 

18 

GR0W7 

-4.7787811814712E+07 

174 

18 

10.34 

40.82 

19 

BANDH 

-1 . 5862801845012E+02 

456 

25 

5.48 

98.82 

20 

E226 

-1.8751929066371E+01 

494 

108 

21.86 

74.55 

21 

STANDATA 

1.2576995000000E+03 

106 

46 

43.40 

18.80 

22 

SCSD1 

8 . 6666666743334E+00 

508 

321 

63.19 

46.96 

23 

GFRDPNC 

6 . 9022359995488E+06 

687 

280 

40.76 

193.26 

24 

BEACONFD 

3 . 3592485807200E+04 

116 

21 

18.10 

13.60 

25 

STAIR 

-2.5126695119296E+02 

415 

42 

10.12 

202.79 

26 

SCRS8 

9. 0429998618888E+02 

609 

171 

28.08 

167.78 

27 

SEBA 

1.571 1600000000E+04 

463 

52 

11.23 

118.27 

28 

SHELL 

1 . 2088253460000E+09 

304 

54 

17.76 

77.81 

29 

PIL0T4 

-2 . 58113926 16949E+03 

1870 

182 

9.73 

827.58 

30 

SCFXM2 

3 . 6660261564999E+04 

858 

130 

15.15 

298.76 

31 

SCSD6 

5 . 0500000078262E+01 

1425 

586 

41.12 

215. 19 

32 

GR0W15 

-1 . 0687094 129358E+08 

435 

50 

11.49 

179.69 

33 

SHIP04S 

1 . 7987 147  004454E+06 

151 

27 

17.88 

36.40 

34 

FFFFF800 

5 . 5567958085232E+05 

1073 

353 

32.90 

341.59 

35 

GANGES 

-1 . 09585 98988428E+05 

703 

234 

33.29 

384.19 

36 

SCFXM3 

5 . 490 1254549751E+04 

1318 

188 

14.26 

692.73 

37 

SCTAP2 

1 . 724807 142857 1E+03 

750 

392 

52.27 

376.31 

38 

GR0W22 

- 1 . 6083433648256E+08 

638 

72 

11.29 

397.37 

39 

SHIP04L 

1 . 7933245379704E+06 

277 

47 

16.97 

67.91 

40 

PILOTWE 

-2 . 720 1032443839E+06 

4982 

572 

11.48 

3520.13 

41 

SIERRA 

1 . 5394390923795E+07 

1317 

542 

41.15 

701.36 

42 

SHIP08S 

1 . 9200982105346E+06 

269 

61 

22.68 

117.89 

43 

SCTAP3 

1 . 4240000000000E+03 

1096 

624 

56.93 

736.06 

44 

SHIP12S 

1 . 4892361344061E+06 

431 

75 

17.40 

259.39 

45 

25FV47 

5 . 5018458882864E+03 

7514 

440 

5.86 

5058.45 

46 

SCSD8 

9.04O9999992546E+02 

4281 

1693 

39.55 

1722.09 

47 

NESM 

1 . 4076079386175E+07 

3067 

0 

0.00 

1338.27 

48 

CZPROB 

2. 1851966988566E+06 

1497 

62 

4.14 

839.22 

49 

PILOTJA 

-6. 1 130540560627E+03 

7398 

547 

7.39 

6466.02 

50 

SH1P08L 

1. 90905521 13891F.+06 

463 

69 

14.90 

236.16 

51 

SHIP12L 

1.4701879193293E+06 

852 

164 

19.25 

589.06 

52 

80BA03B 

9 . 8722736149636E+05 

10769 

921 

8.55 

12652.68 

53 

PILOT 

-5 . 5746058728842E+02 

18441 

4005 

21.72 

75072.32 

Table  2:  Results  from  MINOS  5.3  using  simplified  EXPAND  procedure 
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Problem 

Objective 

Total 

Degen 

Percent 

Solve  t: 

function 

itns 

steps 

degen 

VAX  II  si 

1 

AF1R0 

-4 . 64753142857 14E+02 

6 

3 

50.00 

0.49 

2 

ADLITTLE 

2. 25494963 16238E+05 

92 

14 

15.22 

5.07 

3 

SC205 

-5. 22020612 11707E+01 

117 

13 

11.11 

15.14 

4 

SCAGR7 

-2 . 33 13897523795E+06 

86 

10 

11.63 

7.32 

5 

SHARE 2B 

-4. 1573224074142E+02 

104 

33 

31.73 

7.80 

6 

RECIPE 

-2.6661600000000E+02 

33 

3 

9.09 

2.20 

7 

VTPBASE 

1 . 2983146246136E+05 

59 

17 

28.81 

6.72 

8 

SHARE1B 

-7.6589318579186E+04 

269 

3 

1.12 

25.28 

9 

B0RE3D 

1 . 3730803942085E+03 

159 

62 

38.99 

23.82 

10 

SCORPION 

1.878 1248227381E+03 

103 

38 

36.89 

19.87 

11 

CAPRI 

2.6900129137682E+03 

228 

33 

14.47 

32.19 

12 

SCAGR25 

-1 . 4753433060769E+07 

361 

74 

20.50 

91.79 

13 

SCTAP1 

1 . 4122500000000E+03 

242 

65 

26.86 

37.33 

14 

BRANDY 

1 . 5185098964881E+03 

462 

57 

12.34 

78.95 

IS 

ISRAEL 

-8 . 9664482 186305E+05 

261 

17 

6.51 

38.20 

16 

ETAMACRO 

-7 . 557 1522106445E+02 

501 

115 

22.95 

106.96 

17 

SCFXM1 

1 . 84167590 28349E+04 

386 

66 

17.10 

72.68 

18 

GR0V7 

-4 . 778781 18 147 12E+07 

184 

32 

17.39 

42.67 

19 

BANDM 

-1 . 58628018450 12E+02 

487 

47 

9.65 

107.71 

20 

E226 

- 1.875 192906637 1E+01 

462 

102 

22.08 

72.75 

21 

STANDATA 

1 . 2576995000000E+03 

97 

52 

53.61 

17.47 

22 

SCSD1 

8.6666666743334E+00 

427 

225 

52.69 

38.28 

23 

GFRDPNC 

6 . 9022359995488E+06 

717 

337 

47.00 

206.55 

24 

BEACONFD 

3.3592485807200E+04 

116 

22 

18.97 

14.10 

25 

STAIR 

-2.5126695119296E+02 

364 

40 

10.99 

190.08 

26 

SCRS8 

9. 04299986 18888E+02 

625 

130 

20.80 

177.86 

27 

SEBA 

1.571 1600000000E+04 

417 

44 

10.55 

106.56 

28 

SHELL 

1 . 2088253460000E+09 

310 

54 

17.42 

78.57 

29 

PIL0T4 

-2. 5811392623740E+03 

1452 

141 

9.71 

656.83 

30 

SCFXM2 

3.6660261564999E+04 

880 

138 

15.68 

319.19 

31 

SCSD6 

5 . 0500000078262E+01 

1099 

503 

45.77 

164.71 

32 

GR0V15 

- 1 . 06870941 29358E+08 

446 

65 

14.57 

194.65 

33 

SHIP04S 

1.7987147004454E+06 

149 

25 

16.78 

35.20 

34 

FFFFF800 

5. 5567961145338E+05 

866 

293 

33.83 

281.97 

35 

GANGES 

-1.0958635746320E+05 

679 

218 

32.11 

372.73 

36 

SCFXM3 

5.4901254549751E+04 

1184 

180 

15.20 

632.04 

37 

SCTAP2 

1 . 7248071428571E+03 

680 

386 

56.76 

342.76 

38 

GR0W22 

-1 . 6083433648256E+08 

643 

83 

12.91 

403.74 

39 

SHIP04L 

1.7933245379704E+06 

266 

38 

14.29 

67.03 

40 

PILOTVE 

-2.7201043693969E+06 

5267 

598 

11.35 

3850.05 

41 

SIERRA 

1.5394362183632E+07 

1266 

568 

44.87 

700.02 

42 

SHIP08S 

1. 9200982 105346E+06 

254 

59 

23.23 

113.50 

43 

SCTAP3 

1 . 4240000000000E+03 

840 

502 

59.76 

570.60 

44 

SHIP12S 

1.4892361344061E+06 

445 

87 

19.55 

274.72 

45 

25FV47 

5 . 60 18467790998E+03 

8136 

837 

10.29 

5722.41 

46 

SCSD8 

9 . 0499999992546E+02 

2857 

1251 

43.79 

1174.23 

47 

NESH 

1 . 4076087003981E+07 

2853 

34 

1.19 

1296.87 

48 

CZPROB 

2. 1851966988566E+06 

1503 

130 

8.65 

836.44 

49 

PILOTJA 

-6.11311S2948481E+03 

6141 

422 

6.87 

5496.13 

50 

SHIP08L 

1. 90905521 13891E+06 

467 

67 

14.35 

244 . 25 

51 

SHIP12L 

1.470 1879 193293E+06 

891 

222 

24.92 

621.37 

52 

80BAU3B 

9 . 87231 97704930E +05 

9693 

1068 

11.02 

11768.52 

53 

PILOT 

-5. 5740387782461E+02 

17716 

1624 

9.17 

74443.58 

Table  3:  Results  from  MINOS  5.3  using  maximum-pivot  EXPAND  procedure 
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An  anti-cycling  procedure 


Figure  4:  Comparison  of  iterations  required  by  the  simplex  method  with  different 
steplength  procedures.  The  ratios  ti/1'3  (...)  and  t'2/1'3  ( — )  are  plotted  for  53  test 
problems,  where  (tj,  i2,  t3)  are  the  iterations  for  the  (textbook,  simplified  EXPAND, 
maximum-pivot  EXPAND)  procedures  respectively. 


Figure  5:  Similar  comparison  of  times  required  by  the  simplex  method  with  different 
steplength  procedures. 
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10.  Conclusions 

The  linear  program  (1)  involves  general  constraints  Ax  =  b  and  bounds  l  <  x  <  u. 
The  simplex  method  aims  to  satisfy  Ax  -  b  to  machine  precision,  while  working 
towards  a  solution  that  satisfies  the  bounds  to  some  looser  tolerance  6.  (The  opposite 
is  true  for  certain  other  iterative  methods.) 

The  EXPAND  procedure  was  developed  in  response  to  sporadic  failures  that 
occurred  during  Lustig’s  experiments  with  MINOS  5.1  on  the  same  53  test  problems 
used  here  [Lus87].  We  have  not  experienced  any  failures  since. 

Perhaps  the  main  advance  has  been  in  treating  the  infeasible  blocking  variables 
generated  by  a  Harris-type  ratio  test.  By  retaining  the  infeasible  values  when  such 
variables  become  nonbasic,  we  satisfy  Ax  =  b  to  machine  precision  throughout. 
Only  then  can  we  take  correct  advantage  of  satisfying  bounds  loosely  in  the  manner 
pioneered  by  Harris.  An  important  benefit  is  that  there  is  virtually  no  reversion  to 
Phase  1  after  refactorization — a  common  occurrence  previously  on  ill-conditioned 
problems. 

The  precaution  of  expanding  S  every  iteration  provides  added  theoretical  as¬ 
surance  of  convergence  (given  the  consequent  similarity  to  Wolfe’s  anti-degeneracy 
procedure),  as  well  as  added  practical  assurance  in  the  presence  of  rounding  error. 
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